/*
 * Solver.cpp
 *
 *  Created on: 21.08.2010
 *      Author: jagiella
 */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#include "Solver.h"
#include "Mathematix.h"


// ************ BiCGSTAB ************ //



// x = x - (omega*w + alpha*p)
void calculate_x( float *x, float omega, float *w, float alpha, float *p, int N)
{
	for( int i=0; i<N; i++)
		x[i] -=	(omega*w[i] + alpha*p[i]);
}

// r = r - (omega*v + alpha*s )
void calculate_r( float *r, float omega, float *v, float alpha, float *s, int N)
{
	for( int i=0; i<N; i++)
		r[i] -=	(omega*v[i] + alpha*s[i]);
}

// p = r - beta*(p - omega*s)
void calculate_p( float *x, float *y, float alpha, float *a, float beta, float *b, int N)
{
	for( int i=0; i<N; i++)
		x[i] = y[i] + alpha*a[i] + beta*b[i];
}

// p = r - beta*(p - omega*s)
void calculate_p1( float *x, float *y, float alpha, float *a, int N)
{
	for( int i=0; i<N; i++)
		x[i] = y[i] + alpha*a[i];
}




// x = x - (omega*w + alpha*p)
void calculate_x( double *x, double omega, double *w, double alpha, double *p, int N)
{
#pragma omp parallel for
	for( int i=0; i<N; i++)
		x[i] -=	(omega*w[i] + alpha*p[i]);
}

// r = r - (omega*v + alpha*s )
void calculate_r( double *r, double omega, double *v, double alpha, double *s, int N)
{
#pragma omp parallel for
	for( int i=0; i<N; i++)
		r[i] -=	(omega*v[i] + alpha*s[i]);
}

// p = r - beta*(p - omega*s)
void calculate_p( double *x, double *y, double alpha, double *a, double beta, double *b, int N)
{
#pragma omp parallel for
	for( int i=0; i<N; i++)
		x[i] = y[i] + alpha*a[i] + beta*b[i];
}

// p = r - beta*(p - omega*s)
void calculate_p1( double *x, double *y, double alpha, double *a, int N)
{
#pragma omp parallel for
	for( int i=0; i<N; i++)
		x[i] = y[i] + alpha*a[i];
}

